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Abstract 

Many phenomena, both natural and human-influenced, give rise to signals whose statistical 
properties change under time translation, i.e., are nonstationary. For some practical purposes, 
a nonstationary time scries can be seen as a concatenation of stationary segments. However, the 
exact segmentation of a nonstationary time series is a hard computational problem which cannot be 
solved exactly by existing methods. For this reason, heuristic methods have been proposed. Using 
one such method, it has been reported that for several cases of interest — e.g., heart beat data and 
Internet traffic fluctuations — the distribution of durations of these stationary segments decays with 
a power law tail. A potential technical difficulty that has not been thoroughly investigated is that a 
nonstationary time series with a (scale-free) power law distribution of stationary segments is harder 
to segment than other nonstationary time series because of the wider range of possible segment 
sizes. Here, we investigate the validity of a heuristic segmentation algorithm recently proposed by 
Bernaola-Galvan et al. [Phys. Rev. Lett. 87, 168105 (2001)] by systematically analyzing surrogate 
time series with different statistical properties. We find that if a given nonstationary time series 
has stationary periods whose size is distributed as a power law, the algorithm can split the time 
series into a set of stationary segments with the correct statistical properties. We also find that 
the estimated power law exponent of the distribution of stationary-segment sizes is affected by (i) 
the minimum segment size, and (ii) the ratio R = cr^/ax where is the standard deviation of the 
mean values of the segments, and is the standard deviation of the fluctuations within a segment. 
Furthermore, we determine that the performance of the algorithm is generally not affected by 
uncorrelated noise spikes or by weak long-range temporal correlations of the fluctuations within 
segments. 

PACS numbers: 05.40.-a 
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I. INTRODUCTION 



A stationary time series has statistical properties that do not change under time transla- 
tion ^|. Interestingly, the time series that arise in a large number of phenomena in a broad 
range of areas — including physiologic systems, economic systems, vehicle traffic systems. 



and the Internet 
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ll| — are nonstationary. Thus nonstationarity is a 



property common to both natural and human-influenced phenomena. For this reason, the 
statistical characterization of the nonstationarities in real-world time series is an important 
topic in many fields of research and numerous methods of characterizing nonstationary time 
series have been proposed [l^ . 

One useful approach to quantifying a nonstationary time series is to view it as consisting 
of a number of time segments that are themselves stationary. The statistical properties 
of the segments (i) can help us better understand the overall nonstationarity of the time 
series and (ii) yield practical applications. For example, developing control algorithms for 
Internet traffic will be easier if we understand the statistical properties of these stationary 
segments — a behavior that directly corresponds to the "coarser flow rate" of the current 
network traffic. 

In general, it is impossible to obtain the exact segmentation of a nonstationary time series 

because of the complexity of the calculation. An exact segmentation algorithm requires a 

computation time that scales as 0{N^), where N is the number of points in the time series. 

Hence, such an algorithm is not practical. For this reason, the segmentation of a real-world 

time series must accomplish a trade-off between the complexity of the calculation and the 

desired precision of the result. 

_Bernaola-Galvan and co-workers recently proposed a heuristic segmentation algorithm 

designed to characterize the stationary durations of heart beat time series. In this 

algorithm , the calculation cost is reduced by iteratively attempting to segment the time 

series into only two segments. A decision to cut the times series is made by evaluating a 

modified Student's t-test for the data in the two segments. 

The application of this segmentation algorithm reveals that the distribution of the sta- 
ll 

tionary durations in heart beat time series decays as a power law 13]. Intriguingly, a recent 
analysis of Internet traffic uncovered that the distribution of stationary durations in the 
fluctuation of the traffic flow density also follows a power law dependence Q]. Because 
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these signals have their origin in such diverse contexts, these findings suggest that the power 
law decay of the distribution of the stationary period may be a common occurrence for 
complex time series. Thus, the correct implementation and interpretation of the results 
obtained by the segmentation algorithm is essential in understanding the dynamics of the 
system. In fact, there are many implementation issues concerning the segmentation algo- 
rithm of Ref. [13[ that have not yet been addressed explicitly in the literature, especially 
those concerning the proper estimation of the value of the power-law tail's exponent in the 
cumulative distribution of stationary durations. 

In this paper we systematically analyze different types of surrogate time series to deter- 
mine the scope of validity of the segmentation algorithm of Bernaola-Galvan and co-workers 
jisj ]. In Section II, we briefly explain the segmentation algorithm. In Section III, we present 
results for the dependence of the exponent of the power law tail on the minimum size of 
the segments in the distribution of the stationary durations. In Section IV, we consider the 
effect of the amplitude of the noise and the presence of spike-type noise. In Section V we 
consider long-ranged temporally-correlated noise. Finally, in Section VI, we summarize our 
findings. 



II. IMPLEMENTING THE SEGMENTATION ALGORITHM 



A. The algorithm 

To divide a nonstationary time series into stationary segments we move a sliding 
pointer from left to right along the time series and, at each position of the pointer, compute 
the mean of the subset of the signal to the left of the pointer /xieft and to the right /iright- 
For two samples of Gaussian distributed random numbers, the statistical significance of the 
difference between the means of the two samples, /xi and /xi, is given by the Student's t-test 
statistic {3| 



t 



s 



D 



(1) 



where 



is the pooled variance 'lf)|], si, S2 are the standard deviations of the two samples, and A^i 
and are the number of points in the two samples. 
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Moving the pointer along our time series from left to right, we calculate t as a function of 
the position in the time series. We use the statistic t to quantify the difference between the 
means of the left-side and right-side time series. Larger t means that the values of the mean 
of both time series are more likely to be significantly different, making point tmax? with the 
largest value of t, a good candidate cut point. 

We then calculate the statistical significance P(tmax)- Note that P(tmax) is not the stan- 
dard Student's t-test because we are not comparing independent samples. P(tmax) is numer- 
ically approximated as 



where 7 = 4.19 In — 11.54 and 6 = 0.40 are obtained from Monte Carlo simulations |l3l |. 

is the size of the time series to be split, u = N — 2, and /^(a, b) is the incomplete beta 
function. 

If the difference in mean is not statistically significant — i.e., if is smaller than a threshold 
(typically set to 0.95) — then the time series is not cut. If the difference in means between 
the left and right part of the time series is statistically significant, then the time series is cut 
into two segments as long as the means of the two new segments are significantly different 
from the means of the adjacent segments. If the time series is cut, we continue iterating the 
above procedure recursively on each segment until the obtained significance value is smaller 
than the threshold, or the size of the obtained segments is smaller than a minimum size 
io. We will see that the value of io is one of the parameters controlling the accuracy of the 
algorithm. 

B. Surrogate Time Series 

To investigate the validity of the algorithm, we analyze surrogate time series x{t) gen- 
erated by linking segments with different means. As described in Section I, the cumulative 
distribution of the stationary durations for some real-world time series is characterized by 
a power law decay in the tail, so the probability of finding a segment of size larger than m, 
i.e., the cumulative distribution of segment sizes in our time series, is 




(3) 
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FIG. 1: Surrogate time series constructed according to the procedure described in II. B with 
parameters mo = 20, 7 = 1.0, and R = 1. The time series is drawn in light gray, and the 
stationary segments are represented by a dark gray line. The black line displays the output of 
the segmentation algorithm for £0 = 50. It is visually apparent the black line provides a adequate 
coarse-grained description of the surrogate time series. Note also that the algorithm cannot extract 
the smallest segments because of the restrictions on resolution for m < Iq. 

where mo is the minimum size of a segment. 

We generate time series with a power law distribution of segment sizes by the following 
procedure: 

1. draw from the interval [mo, +00] a sequence of segment lengths {nii} with distribution 
given by Eq. (jH); 

2. draw from the interval [0,1] a sequence of mean time series values Xi with uniform 
probabihty; 

3. draw from the interval [-v^o'e,\/3o'e] a sequence of fluctuation values ei{ki), for fcj = 
1, ■ ■ ■ , mj, with uniform probability. 



The resulting time series is given by 



x(t) = Xi + tiiki) 



where i is such that 



i-l 



< t < ^m^-, 
j=o i=0 



(5) 



(6) 



and ki is such that 

i-l 
j=0 

To quantify the level of the noise, we define the ratio 

R^-, (8) 

where ax is the standard deviation of the mean of the segments and is the standard 
deviation of the fluctuations within a segment. For x uniformly distributed in the interval 
[0,1], a, = 0.3. 

For each set of parameters {R, 7, mo) we generate ten time series, each with 50,000 data 
points. Note that knowing a priori the value of mo in a real-world time series is unlikely, 
but in order to test the algorithm in a consistent way, we consider in the following section 
mo > 20 because that is the resolution limit for the algorithm (see also Appendix A). 

Figure HI displays a surrogate time series, the corresponding stationary durations, and the 
result of the segmentation algorithm with = 50 and mo = 20. The segments obtained by 
the segmentation algorithm do not exactly match the stationary segments in the surrogate 
time series but the figure strongly suggests that the algorithm provides the correct coarse- 
grained description of the time series. As expected, the segmentation algorithm cannot 
extract segments with size m < Cq. 



III. ACCURACY OF THE SEGMENTATION ALGORITHM 



A. Dependence on £0 and mo 

Figure El^a) displays the cumulative distribution of segment sizes, which is the probability 
of finding a segment with size larger than £ for a surrogate time series split for different 
values £o- The cumulative distributions of the size of the stationary segments cut by the 
segmentation algorithm for surrogate time series are well fit by a power law decay. For all 
cases, we find 

P£o,mo(>^)~^"^^'"'"°^ (9) 

with the same exponent value 7 1.0 indicating that the segmentation algorithm 

splits the nonstationary time series into segments with the correct asymptotic statistical 
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FIG. 2: (a) Cumulative distribution P{£) of segment sizes larger than £ for surrogate time series 
with 50,000 data points; 7 = 1.0, R = 1 and mo = 20. The dotted line indicates the input segment 
size distribution for the surrogate time series. The slope of the hne is 1.0 for 20 < m < 4,000. 
The other curves show P{i) for different values of Iq. For Iq = 20, the curve is well described in 
the range 100 < £ < 4,000 by a power law with 7 « 1.0. The distribution does not decay as a 
power law for i < 100 due to the fact that the segmentation algorithm cannot split a time series 
with a number of points insufficient to perform the Student's f-test. For Iq = 400, P{i) decays 
as a power law for 1,000 < £ < 8,000. The segmentation algorithm correctly splits all segments 
in the surrogate scries with size £ > 1,000. (b) Dependence of 7 on £q and mg. The mean and 
the standard deviation of the exponent for the original time series is 1.05 ± 0.24, shown by the 
black solid and dotted lines. The error bars show the standard deviation of the estimates 7. For 
£o/mo < 5, we find 7 1. Thus, the algorithm accurately estimates exponents in this region. 
However, the values of the exponent are close to 1.3 for io/mo > 10, meaning that mo <C £0 leads 
to an overestimation of 7. 



properties. However, the range of scales for which we observe a power law decay with 7 pa 7 
depends strongly on the selection of £9. 

For i greater than about 5£o, the tails of the distributions are close to power law decays 
with 7 fa 7 = 1.0. Moreover, all P{£) follow the original curve for £ > 1,000, i.e., the 
algorithm correctly identifies the large segments independently of the selection of £0. For 
£ < 5£o, the distributions do not follow power law decays. The origin of this behavior hes in 
the fact that (i) for io = 20 = 0(mo), there are not enough data points to reliably perform 
the Student's t-test, so one cannot reasonably expect any statistical procedure to be able to 
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TABLE I: Estimated exponent 7, as defined by Eq. (5) for 7 = 1.0. The mean and the standard 
deviation of the exponents are calculated for the ranges indicated inside parenthesis using Eq. HB1() . 
The column labeled "input" presents exponent estimates obtained from the segment lengths used 
to generate the surrogate time series. We find 7 ~ 7 = 1.0 for niQ < 5 £q. 



\ ^0 

mo \ 


10 


20 


50 


100 


400 


input 


20 


1.1 ±0.3 

{£ > 20) 


1.1 ±0.3 

{£ > 20) 


1.2 ±0.3 

{£ > 50) 


1.3 ±0.4 
{£ > 110) 


1.3 ±0.3 
{I > 1000) 


1.1 ±0.3 


50 


1.0 ±0.3 

{£ > 50) 


1.0 ±0.3 

{£ > 50) 


1.0 ±0.3 

{i > 100) 


1.1 ±0.3 

{£ > 110) 


1.1 ±0.5 

(l > 1000) 


1.0 ±0.2 


100 


1.0 ±0.3 

{£ > 100) 


1.0 ±0.3 
{i > 100) 


1.0 ±0.3 

{£ > 100) 


1.0 ±0.3 

{£ > 100) 


1.2 ±0.3 
(i > 1000) 


1.1 ±0.3 



extract those short segments, and (ii) for ip ^ mo, one fails to extract many of the stationary 
segments. The reason for the latter is that the value of io is in this case considerably larger 
than the size of the shortest segments in the time series, so the algorithm is forced to merge 
a number of short segments into longer ones with size greater than £q. This process gives 
rise to an excess of segments with sizes between ip and 2£o- 

Table |l] shows the mean and standard deviation of the estimated exponent value 7 calcu- 
lated from surrogate time series for several values of mo and io (see Appendix A for details on 
how to to estimate 7). Our results indicate that 7 depends on both mo and Eq: If io ^ mo, 
7 overestimates 7, while if io ~ 0(mo), the algorithm correctly estimates the value of the 
exponent 7; cf. Fig. |H^b). 

B. Dependence on 7 

Next, we focus on the dependency of the accuracy of the segmentation algorithm on the 
value of 7. Figure a) displays the cumulative distribution of segment sizes for surrogate 
time series generated with 7 = 2.0. A challenge for the segmentation algorithm is that 
Eq. (j3]) indicates that the probability of finding segments with size shorter than 4mo is 90% 
for 7 = 2.0, which can lead to the "aggregation" of several consecutive segments of small 
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FIG. 3: (a) Cumulative distribution of segment sizes for 7 = 2.0, i? = 1 and niQ = 20. For £q = 20 
and 100 < £ < 1,000, the exponent 7 of the power law is close to 7. For io = 100, we also 7 7. 
However, for = 400, the algorithm fails to split the time series correctly for i < 1,000. Moreover, 
note that even though the exponent estimate is correct, the algorithm yields segments that are 
longer than the ones in the surrogate time series, (b) Dependence of 7 on Iq and tuq for 7 = 2.0. 
For io/rriQ < 4, the algorithm yields segments with the correct statistical properties. 

size into a single longer segment. 

As we found for 7 = 1, the tails of the distributions follow power law decays for large 
i, showing that the algorithm yields segments with the proper statistical properties. In 
Fig. ini^b), we show the dependence of 7 on io and mo for 7 = 2.0. We find that for 
io/rriQ < 4 the algorithm extracts segments with distributions of lengths that decay in the 
tail as power laws with exponents that are quite close to 2.0 . 

In Tables ITTl and ITTTl we report the values of 7 for 7 = 2.0 and 7 = 3.0, respectively. We 
find that for small mo and large 7, one over-estimates 7. Thus, we surmise that for 7 > 3.0, 
it becomes impractical to estimate 7 accurately, except for extremely long time series. This 
fact is not as serious a limitation as one may think because for large 7 it is always difficult 
to judge whether a distribution decays in the tail as an exponential or as a power law with 
a large exponent. 
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TABLE II: Estimated exponents 7 for 7 = 2.0. The mean and the standard deviation of the 
exponents are calculated for the ranges indicated inside parenthesis using Eq. HB1|) . The column 
labeled "input" presents exponent estimates obtained from the segment lengths used to generate 
the surrogate time series. We find 7 7 = 2.0 for mo < 4 ip. 



mo \^ 


10 


20 


50 


100 


400 


input 


20 


2.2 ±0.6 

{£ > 100) 


2.2 ±0.5 
{i > 100) 


2.2 ±0.5 

{£ > 100) 


2.4 ±0.7 

{£ > 200) 


3.8 ±0.7 

{£ > 800) 


2.2 ±0.1 


50 


2.2 ±0.6 

{£ > 100) 


2.2 ±0.6 
{i > 100) 


2.2 ±0.6 

{£ > 100) 


2.5 ±0.6 

{£ > 200) 


3.8 ±2.1 
{£ > 1000) 


2.1 ±0.4 


100 


1.9 ±0.5 

{£ > 100) 


1.9 ±0.6 
{i > 100) 


1.9 ±0.6 

{£ > 100) 


2.0 ±0.5 

{£ > 200) 


2.4 ±0.4 
{£ > 1000) 


2.0 ±0.4 



TABLE III: Estimated exponents 7 for 7 = 3.0. The mean and the standard deviation of the 
exponents are calculated for the ranges indicated inside parenthesis using Eq. HB1|) . The column 
labeled "input" presents exponent estimates obtained from the segment lengths used to generate 
the surrogate time series. For this value of 7, one finds that a small mo leads to a clear under- 
estimation of the 7. Note that the standard deviation of 7 becomes larger, indicating the difficulty 
in obtaining 7 accurately. Also noteworthy is the fact that because 7 is so large, the range of 
segment lengths m drawn becomes much reduced. This implies that if one sets £q = 400 one is 
unable to properly estimate 7. 



\ ^0 
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10 


20 


50 
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400 


input 


20 


2.4 ±0.5 


2.3 ±0.5 


2.5 ±0.5 


2.7 ±0.5 


4.3 ± 1.4 


3.1 ±0.5 




{£ > 100) 


{£ > 100) 


{£ > 100) 


{£ > 200) 


{£ > 1000) 




50 


3.0 ±0.8 


3.0 ±0.8 


3.0 ±0.9 


3.4 ±0.6 


5.9 ± 1.3 


3.0 ±0.6 




{£ > 100) 


{£ > 100) 


{£ > 100) 


{£ > 200) 


{£ > 1000) 




100 


2.8 ±0.9 


2.8 ±0.8 


2.8 ±0.9 


2.8 ±0.8 


5.2 ± 1.7 


2.9 ±0.7 




{£ > 200) 


{£ > 200) 


{£ > 200) 


{£ > 200) 


{£ > 1000) 
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IV. ROBUSTNESS OF THE ALGORITHM WITH REGARD TO NOISE 



A. Amplitude of fluctuations around a segment's mean 



Another factor that may affect the performance of the segmentation algorithm of 
Bernaola-Galvan and co-workers is the amphtude of the fluctuations within a segment. 
It is plausible that greater noise amplitudes will increase the difficulty in identifying the 
boundaries of neighboring segments. Thus, we next analyze the effect of the amplitude of 
the noise for surrogate time series. 

Figure E] demonstrates that for large R, the segmentation algorithm yields few short 
segments. This result arises from the concatenation of neighboring segments with means 
that become statistically indistinguishable due to the large value of a^. 

3 1 1 1 1 1 1 1 1 1 1 — 
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FIG. 4: Surrogate time series for different amplitudes of the fluctuations within the segments. The 
surrogate time series were generated with mo = 20, 1^ = 20, 7 = 1.0, and (a) R = 1.0, (b) R = 2.0, 
(c) R = 3.0, and (d) R = 4.0. It is visually apparent that the segmentation algorithm yields longer 
segments as R increases. This fact arises from the fact that the statistical test cannot distinguish 
two neighboring segments with close — relative to — means. 
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FIG. 5: Cumulative distribution of segment size obtained with the segmentation algorithm for 
different amplitudes of the noise. As in Fig.^ the time series analyzed were generated for rriQ = 20, 
7 = 1.0, and (a) R = 1.0, (b) R = 2.0, (c) R = 3.0, and (d) R = 4.0. For < 3, the tails of 
the distributions decay as power laws. For R = 4.0, it is difficult to discriminate whether our the 
tails of distribution conform to exponential or power law decays. Note that as R increases, the 
dependence of the functional form of the distributions on decreases appreciable. 



We show in Fig. El the cumulative distributions of segment sizes for 7 = 1.0, mo = 20, 
and for different values of R. For large ^ and i? < 3, we find 7^7, while for i? = 4 the 
algorithm becomes ineffective at extracting the stationary segments in the time series. It is 
visually apparent that for i? = 4.0 the fluctuations within a segment are so much larger than 
the jumps between the means of the stationary segments that the segmentation becomes 
totally unable to parse the different segments; cf. FigHjl. 

In Fig. ini we show 7 for different values of R. For 7 = 1.0, we estimate 7 7 for 
< i? < 4. In contrast, for 7 = 2.0, 7^7 only for < i? < 1.5. When R increases, the 
segmentation algorithm is unable to cut the segments because the greater amplitude of the 
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FIG. 6: Dependence of 7 on the ratio R = o^jox for mo = 20, 1^ = 20 and three distinct values 
of 7: (a) 7 = 1.0, (b) 7 = 2.0, and (c) 7 = 3.0. The different curves in each plot correspond to 
different values of a^. For 7 = 1.0, we find 7 = 1.2 for < R < 4.0, indicating that the algorithm 
is robust against increases in R. For 7 > 1.0, we find that the impact of an increasing R on the 
performance of the algorithm becomes more and more marked. Specifically, for 7 = 2.0, we find 
7 ~ 7 for R < 1.5, while for 7 = 3.0, we find 7 7 for R < 0.6. 

fluctuations inside a segment decreases the significance of the differences between regions of 
the time series. This effect yields very large segments, which results in very small estimates 
of 7. This effect is even stronger for 7 = 3.0, for which we find 7^7 only for < i? < 0.6. 

B. Spike noise 

Next we analyze the effect of spike noise on the performance of the segmentation algo- 
rithm. We generate surrogate time series as before and then for each t make, with probability 
p, x{t) = 2. The effect of this procedure is illustrated in Fig|7|for four distinct values of 
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FIG. 7: Surrogate time series with uncorrelated spike noise for niQ = 20, 7^ = 0.3, Iq = 20, 7 = 1.0 
and (a) p = 0, (b) p = 0.01, (c) p = 0.05, and (d) p = 0.1. For all values of p considered, the 
application of the segmentation algorithm yields segments that, in a coarse-grained way, match 
well the segments in the surrogate data the. 

p. The figure also suggests that the segmentation algorithm yields a good coarse-grained 
description of the surrogate time series for p as large as 0.1. This result suggests that the 
algorithm is robust to the existence of uncorrelated spike noise in the data (Fig. 



V. CORRELATED NOISE 

In this Section, we investigate the effect of long-range correlations in the fluctuations 
around the segment's mean on the performance of the segmentation algorithm. This study 
is particularly important because real-world time series often display long-range power law 
decaying correlations. 
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FIG. 8: Cumulative distribution of segment sizes for surrogate time series with 7 = 1.0, = 0.3, 
mo = 20, and for different densities p of the spike noise: (a) p = 0.0, (b) p = 0.01, (c) p = 0.05, 
and (d) p = 0.1. The value of x{t) is set at 2.0 for a spike. It is visually apparent that even for 
p = 0.1 all the distributions decay as power laws for i > 200 and that the slopes in the log-log plots 
are similar to the slope of the input distributions, i.e., the segmentation algorithm yields segments 
with the correct statistical properties even in the presence of strong spike noise. 

A. Segmentation of correlated noise with no segments 

We generate a temporally-correlated noise whose power spectrum decays as S{f) ~ 
jisl ]. The surrogate time series consists of 60,000 points, with mean 0. Figure 01 displays the 
cumulative distribution of segment lengths. The curves show the different noise correlations: 
P = 0.3, 0.5, and 1.0 (1//). Note that we have confirmed only one long segment for (3 = 0.0, 
as one expects. On the other hand, for small (3, we still observe a longer segment, although 
the curves decay rapidly for large (3. Thus, the segmentation algorithm divides a correlated 
noise for large (3 into many small pieces of stationary segments, because the given time series 
is already nonstationary for this case. Also, it is hard to determine the functional form of 
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the plots for all cases. Even if the plots are assumed to follow power laws, the relationship 
between (3 and 7 is quantitatively obscure. 
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FIG. 9: Segmentation of Guassian distributed long-range correlated noise with different exponent 
of the power law in the power spectrum density for the noise; (a) /3 = 0.3, (b) f3 = 0.5, and (c) 
P = 1.0. We generate time series with 60,000 data points according to the modified Fourier filtering 



method (FFM) [1^. It is hard for all the curves to determine the functional form; a power law 
or an exponential distribution. There exists only one long segment for /? = 0. Then, for larger f3, 
the segmentation algorithm splits the time series into many smaller pieces of stationary segments. 
This is because that the time series with larger /3 has a higher probability that a part of the time 
series is nonstationary. However, a quantitative relationship between 7 and /? is not clear for all 



the case of the correlated noise with no segments. 



B. Segmentation of correlated noise with segments 

Here, in order to show the ability of the identification of segments, we analyze surrogate 
time series concatenating segments with long-range correlated random variables whose power 
spectrum is followed by a power law with an exponent /3. We show in Figs. Il()r a)-(c) typical 
surrogate time series generated with 7 = 1.0, mo = 20, R = 1.0, £0 = 20, and different 
temporal correlations: (a) P = 0.1, (b) P = 0.3, and (c) (3 = 0.5. The figure suggests that 
the segmentation algorithm can correctly parse the short segments but that long segments 
get cut multiple times, especially for (3 = 0.5. This result is to be expected because the 
strong correlations in the noise lead to marked changes in the mean. 

In Figs. ITT¥ d)-(f). we show the cumulative distributions of segment sizes for (5 = 
0.1, 0.3, and 0.5. The data confirm quantitatively the visual impression gained from 
Figs, unfa)- (c), i.e., that longer segments get cut multiple times. In particular, for /5 = 0.5, 
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the distributions clearly deviate from the power law, independent of the selection of io- How- 
ever, this fact should not be seen as a shortcoming of the algorithm; for large (3, correlated 
noise for a long segment is already nonstationary, so that the algorithm is indeed cutting a 
nonstationary signal into stationary durations. 
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FIG. 10: Surrogate time series with the long-range correlated noise and the result of the segmen- 
tation for (a) /3 = 0.1, (b) (3 = 0.3, and (c) (3 = 0.5. The time series was concatenated by segments 
with the long-range correlated noise generated by the modified Fourier filtering method (FFM) [l^ 
for 7 = 1.0, rriQ = 20, = 20, and R = I. For /? = 0.1 and 0.3, the results of the segmentation 
algorithm resemble the segments in the surrogate time series. However, for /3 = 0.5, it is visually 
apparent that long segment are cut multiple times, indicating that the algorithm judges the noise 
within a segment as a nonstationary time series. This is because that for large /?, the correlated 
noise in a segment cannot be stationary. The same is not true for short segments, which are iden- 
tified correctly under the above conditions. Cumulative distributions of the segument size for (d) 
f3 = 0.1, (e) P = 0.3, and (f) (3 = 0.5. The distributions confirm the visual impression obtained 
form (a)-(c). In summary, the distributions decay faster than a power law when /3 > 0.3, although 
they still follow a power law dependence for smaller /?. 
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VI. DISCUSSION 



In this paper we analyzed nonstationary surrogate time series with different statistical 
properties in order to investigate the validity of the segmentation algorithm of Bernaola- 

n 

Galvan and co-workers Our results demonstrate that this heuristic segmentation al- 

gorithm can be extremely effective in determining the stationary regions in a time series 
provided that a few conditions are fulfilled: First, one must have enough data points in the 
time series to yield a large number of segment lengths, otherwise one will not be able to 
reach the aymptotic regime of the tail of the distribution of segment sizes. 

Second, the ratio of the amplitude of the fluctuations within a segment to the typical jump 
between the means of the stationary segments must be relatively small (less than about 0.6) 
in order for one to trust the output of the segmentation algorithm. This concern contrast 
with the case of spike noise in the data which affects the performance of the segmentation 
algorithm only weakly. 

Finally, if there are long-range temporal correlations of the fluctuations around the mean 
of the segment, then the segmentation algorithm correctly cuts the time series into the 
stationary segments for small (3. However, for (3 > 0.3, the fluctuations inside long segments 
become nonstationary, which results in the algorithm detecting many "stationary" durations 
inside these long segments. 

Our analysis provides a number of clear guidelines for using effectively the segmentation 
algorithm of Bernaola- Galvan et al. jl^l : 

1. One must perform the segmentation for a number of different values of ^Q in order to 
identify the region for which the tails of the distributions of segment sizes reach the 
asjTiiptotic scaling behavior. (Note: If 7 is large, than the estimation error can be 
quite considerable, especially if uiq is small.) 

2. One must calculate the ratio between the standard deviation of the mean value of 
the segment and the standard deviation of the fluctuations within a segment after 
performing the segmentation. If the R than 0.6, then there is the possibility that 7 is 
considerably under-estimating the true value of 7. 
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FIG. 11: Algorithm's resolution with regard to io and m. We generate surrogate time series are 
constructed by the alternated concatenation of two types segments with fixed size mo: one with 
mean 0.0 and standard dewviation 0.29 and another with mean 1.0 and standard deviations 0.29. 
We plot Q, which is defined in Eq. (|A1|) . as a function of ni for diff'erent values of £q. Our results 
indicate that the segmentation algorithm is not able to extract the stationary segments for m < 20 
or £o < 10. Note also that £o <m for correct segmentation to occur. 



APPENDIX A: PERFORMANCE OF THE ALGORITHM FOR FIXED SEG- 
MENT SIZES 

As a special case, we analyze the time series whose stationary duration is fixed in order 
to discuss the minimum resolution of the segment algorithm. We concatenate segments 
with constant lenght m with alternating means of 0.0 and 1.0. We then add fluctuations 
to those segments with a standard deviation of 0.29. We define the fraction of successfully 
split segments 

^ _ Number of segments correctly identified by the algorithm (Al) 
Number of segments in surrogate data ' 

where Q = 1.0 corresponds to perfect segmentmentation. 
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We plot in Fig. ^2 Q as a function of m for different values of io . For m < 20, the 
segmentation algorithm does not yield the correct segments in the surrogate time series even 
though the segment's means are quite different. This result suggests that the resolution of 
the segmentation algorithm is ~ 20. We also find that for = 5, the algorithm splits 
the time series into too many segments. This result suggests that for optimal performance 
io > 10. 



APPENDIX B: ESTIMATION OF 7 

To estimate the exponent 7 of the power laws, we calculate the power law exponent at a 
segment size in from a small region around in 

,_ log(P(C+.))-log(P(C)) 

log(4+.) - log(4) ' ^ ' 

where we set the region to r = 5. Figure IT^ shows the behavior of the power law exponent 
over i corresponding to Fig. Efa); the black lines are the behavior of the exponent at seg- 
ment size i, and the broken lines are the behavior of the exponent of the power law in the 
distribution of the surrogate time series for comparison. 

The range of in for the calculation is between a starting point (where the exponential 
decay disappears) and 10,000. For io = 20, the value of the exponent decreases for i < 100, 
which corresponds to the exponential decay as shown in Fig. |21^a). However, the value of 
the exponent fluctuates around 1.0 for i > 100, indicating that the algorithm can reproduce 
the correct statistical behavior. For io = 50, the value of the exponent is close to 1.0 for 
50 < i < 10,000. Moreover, for io = 100 and 400, the curves decay quickly for the smaller 
i, and the value of the exponent tends to be overestimated. 
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FIG. 12: The behavior of the exponents ^{in) of the power law in segment size distribution for 
different 4; (a) 4 = 20, (b) Iq = 50, (c) io = 100, and (d) Iq = 400. mo = 20, 7 = 1.0. The 
exponent of the power law is calculated by Eq. HBlj) . and plotted by black lines in the figure. 
The dotted gray lines also indicate the exponent of the surrogate time series for comparison. For 
Iq = 20, the observed curve resembles the original one for i > 100 and stays around 1.0, though 
the curve indicates an exponential decay for i < 100. For Iq = 100 and 400, the segmentation 
results can follow the original ones for larger i; however, the estimated values of the exponent are 
somewhat overestimated for smaller i 
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